Обращение преобразования Миуры

В.Э. Адлер. Классические интегрируемые системы, весенний семестр 2020

В этом файле показано, как можно реализовать символьные вычисления, связанные с оператором полной производной, эволюционным дифференцированием, вариационной производной и алгоритмом интегрирования по частям. Это иллюстрируется на примере задачи о законах сохранения КдФ; здесь мы будем строить их при помощи преобразования Миуры. Определённые в этом файле функции неоднократно пригодятся и в дальнейшем.

Что такое символьные вычисления

Есть несколько уровней работы с дифференциальными уравнениями. При численном построении решения (например, по схеме Забуски-Краскала из первой лекции) обрабатываются массивы вещественных чисел. Тут без компьютера никак не обойтись.

При построении точных решений типа солитонных компьютер тоже помогает – делает вычисления с такими функциями как, например, cosh x и рисует графики.

Есть и более абстрактные вычисления. Пусть задано некоторое ДУ на функцию u(x,t). Можно работать с этой функцией и её производными просто как с некими символами, формальными переменными. При таком подходе мы фактически работаем не с решениями уравнения, а с самим уравнением; всё, что у нас есть – это оно само и его дифференциальные следствия, то есть, соотношения, которые получаются при дифференцировании по x и t (в предположении, что u∈(C^∞)). В теории интегрируемых систем вычисления этого типа встречаются очень часто и играют важную роль. В принципе, тут некоторое время можно обходиться без компьютера. Чтобы как следует понять суть вычислений, которые у нас будут, я советую некоторую их часть делать на бумаге. Например, ниже приведен алгоритм построения законов сохранения для уравнения КдФ; попробуйте несколько первых законов сохранения найти и проверить вручную. Через несколько лекций будут похожие вычисления с высшими симметриями, тоже несколько первых нетрудно найти вручную. Вы увидите, что выкладки простые, но их объем быстро растет в зависимости от порядка по производным. Компьютер просто помогает сделать на несколько шагов больше.

Основные операции

Рассмотрим эволюционные уравнения с одной полевой переменной u=(Subscript[u, 0]):

(1) 	(Subscript[u, t])=a(x,(Subscript[u, 0]),(Subscript[u, 1]),...,(Subscript[u, n])),       (Subscript[u, n]):=(∂_x^n)(u).

Основная идея заключается в том, чтобы считать x,t и все (Subscript[u, n]) независимыми (“динамическими”) переменными. Мы рассматриваем всевозможные (гладкие) функции от конечного числа этих переменных – например, саму правую часть a, правые части высших симметрий, плотности законов сохранения; в общем, всё что угодно. То обстоятельство, что (Subscript[u, n]) имеют смысл производных от u, кодируется в правилах дифференцирования таких функций. Чтобы их задать, достаточно определить следующие функции (не забудьте активировать эту ячейку, для этого нужно поместить на неё курсор и нажать Shift-Enter):

vars[f_]:=Union[Cases[f,_u,{0,∞}]]
dif[f_]:=Plus@@(D[f,#]dd[#]&/@vars[f])
dx[f_]:=D[f,x]+dif[f]/.dd[u[k_]]:>u[k+1]
dx[f_,k_]:=Nest[dx,f,k]
star[g_,f_]:=dif[g]/.dd[u[k_]]:>dx[f,k]
vard[f_]:=Expand[Plus@@((-1)^#[[1]]dx[D[f,#],#[[1]]]&/@vars[f])]

Остальные ячейки в этом разделе можно не запускать, они служат лишь для тестирования и пояснения того, что означают эти определения.

Функция vars[f] возвращает список переменных (Subscript[u, n]) в выражении f (если мы применяем его не упростив f, то этот список может оказаться более широким, чем нужно, это не страшно):

f=(((u[1]+u[13]))^2)-((u[1])^2)-2u[1]u[13]+b[x]
vars[f]
Expand[f]
vars[%]

Функция dif[f] возвращает дифференциал выражения, причем d (Subscript[u, n]) обозначается dd[u[n]]. К самим функциям dif и dd мы обращаться не будем, но они используются в определении функции dx[f], которая возвращает полную производную

(Subscript[D, x])(f)= (∂_x)(f)+(Subscript[u, 1])(∂_0)(f)+... +(Subscript[u, n+1])(∂_n)(f),

где f=f(x,(Subscript[u, 0]),(Subscript[u, 1]),...,(Subscript[u, n])) и (∂_k) обозначает частную производную по (Subscript[u, k]). Обратим внимание, что мы не пытаемся упростить выражение по ходу вычислений. Практика показывает, что это лучше делать в конце, при помощи подходящей команды (Expand, Together, Factor, Collect, Simplify, ...). Если окончательный результат хочется увидеть в индексных обозначениях, это тоже можно сделать (но в процессе вычислений использовать индексы неудобно).

f
dif[f]
dx[f]
Expand[dx[f]]
Expand[dx[f,2]]
%/.u[n_]:>(Subscript[u, n])

Функция star[g,f] возвращает

(SubStar[g])(f)=(Subscript[∇, f])(g)=f(∂_0)(g)+... +(Subscript[D, x]^m)(f)(∂_m)(g) ,

где g=g(x,(Subscript[u, 0]),(Subscript[u, 1]),...,(Subscript[u, m])), то есть, результат дифференцирования выражения g в силу уравнения (Subscript[u, t])=f.

star[u[0],u[3]+6u[0]u[1]]
star[u[1],u[3]+6u[0]u[1]]
star[u[0]u[2],u[3]+6u[0]u[1]]

Функция vard[f] возвращает вариационную производную (= оператор Эйлера) от f:

δ f/δ u=E(f)= (∂_0)(f)-(Subscript[D, x])((∂_1)(f))+(Subscript[D, x]^2)((∂_2)(f))-...+(((-(Subscript[D, x])))^n)((∂_n)(f)),

где f=f(x,(Subscript[u, 0]),(Subscript[u, 1]),...,(Subscript[u, n])) и (∂_k) обозначает частную производную по (Subscript[u, k]).

vard[u[3]+6u[0]u[1]]
vard[u[2]g[u[0]]]

Теорема. Выражение является полной производной по x тогда и только тогда, когда его вариационная производная равна 0:

E(f)=0   ⇔  ∃g:f=(Subscript[D, x])(g).

Напишите здесь любое выражение вместо f (ну, не очень сложное, такое, чтобы его смогла осилить машина!) и проверьте, что в последней строчке получится 0.

f=g[u[3]+u[2]]((u[4])^2)Sin[x]Cos[u[2]]
dx[f]
Expand[vard[%]]

Преобразование Миуры

Рассмотрим модифицированное уравнение КдФ, точнее, (SuperMinus[мКдФ]):

(2)	(Subscript[v, t])=(Subscript[v, x x x])-6((v^2)+λ)(Subscript[v, x]),

с произвольным параметром λ. Определим переменную

(3)	u=-(Subscript[v, x])-(v^2)-λ.

Дифференцируя u по x и по t, с заменой (Subscript[v, t]) из (2), нетрудно проверить (сделайте это вручную), что u удовлетворяет уравнению КдФ

(4)	(Subscript[u, t])=(Subscript[u, x x x])+6u (Subscript[u, x]).

Чтобы проверить это с помощью определенных выше операций, мы временно поменяем ролями u и v (так как правила дифференцирования у нас заданы для u):

ut=u[3]-6(((u[0])^2)+λ)u[1];
v=-u[1]-((u[0])^2)-λ;
-star[v,ut]+dx[v,3]+6v dx[v]
Expand[%]

Замечание. У мКдФ есть две вещественные формы, отличающиеся знаком нелинейного члена:

(Subscript[v, t])=(Subscript[v, x x x])±6((v^2)+λ)(Subscript[v, x]).

Очевидно, они связаны заменой v->i v, λ->-λ, но она не годится, если мы следим за вещественностью решений. С этой точки зрения, свойства этих двух уравнений довольно сильно отличаются (например, у (SuperPlus[мКдФ]) есть решения в виде бризеров, а у (SuperMinus[мКдФ]) нет; представления нулевой кривизны имеют разный вид). Вещественное преобразование Миуры (3) определено только для (SuperMinus[мКдФ]). Однако, сейчас это для нас не важно, и вместо (2), (3) мы могли бы рассматривать уравнения

(Subscript[v, t])=(Subscript[v, x x x])+6((v^2)+λ)(Subscript[v, x]),      u=-i (Subscript[v, x])+(v^2)+λ:

ut=u[3]+6(((u[0])^2)+λ)u[1];
v=-I u[1]+((u[0])^2)+λ;
Expand[-star[v,ut]+dx[v,3]+6v dx[v]]

Преобразование (3) это пример так называемых дифференциальных подстановок. В отличие от точечных замен переменных (таких, как растяжение переменных, переход в подвижную систему координат, переход от декартовых координат к полярным и т.п.) это преобразование необратимо: нет такой формулы, которая выражает v как функцию от конечного числа переменных (Subscript[u, n]). Тем не менее, такие односторонние преобразования тоже могут быть очень полезными. Во-первых, иногда уравнение может быть линеаризовано с их помощью. Вот простой пример, подстановка Коула-Хопфа. Здесь уравнение для u – линейное уравнение теплопроводности, и по любому его решению мы строим решение v нелинейного уравнения Бюргерса  (тоже проверьте на бумаге, это совсем просто):

ut=u[2];
v=u[1]/u[0];
Expand[-star[v,ut]+dx[v,2]+2v dx[v]]

Во-вторых, даже если подстановка не приводит к упрощению уравнения (непонятно, чем (2) лучше чем (4)), её, возможно, удастся как-то использовать для построения частных решений. Действительно, на одной из следующих лекций мы увидим, что преобразования Миуры можно соединить в цепочку, позволяющую строить многосолитонные решения.

В-третьих, прямо сейчас мы увидим, что из формулы (3) можно очень просто получить рекуррентную схему для вычисления законов сохранения уравнения КдФ. Подобный приём срабатывает и во многих других примерах.

Обращение преобразования Миуры

Хотя выше было сказано, что преобразование (3) необратимо, мы всё же можем выразить v в виде ряда по обратным степеням параметра (((-λ))^(1/2)), причём коэффициенты ряда будут многочленами от (Subscript[u, n]). Сходимость ряда нас не интересует – это просто формальный ряд, производящая функция. Положим λ=-(z^2) и перепишем преобразование (3) в виде

(5)	(v^2)+(Subscript[v, x])=(z^2)-u.

Подставим сюда разложение

(6)	v=-z+(Subscript[V, 0])+(Subscript[V, 1])/z+(Subscript[V, 2])/(z^2)+...

(индекс у (Subscript[V, k]) это просто номер, не производная), это даст

(z^2)-2z (Subscript[V, 0])+(-2(Subscript[V, 1])+(Subscript[V, 0]^2))+(-2(Subscript[V, 2])+2(Subscript[V, 0])(Subscript[V, 1]))/z
           + (-2(Subscript[V, 3])+2(Subscript[V, 0])(Subscript[V, 2])+(Subscript[V, 1]^2))/(z^2)+ ...  +(Subscript[V, 0,x])+(Subscript[V, 1,x])/z+(Subscript[V, 2,x])/(z^2) +...     =  (z^2)-u.

Отсюда следует, что (Subscript[V, 0])=0, (Subscript[V, 1])=u/2,  а на остальные коэффициенты получаем рекуррентные соотношения

2(Subscript[V, n+1]) =(Subscript[V, 1])(Subscript[V, n-1])+...+(Subscript[V, n-1])(Subscript[V, 1])+(Subscript[V, n,x])=0,   n=1,2,3,... .

Вычислим несколько первых коэффициентов (скажем, 13, хотя машина может гораздо больше; на бумаге продублируйте, пока не надоест).

M=13;
V[1]=u[0]/2;
Do[V[n+1]=Factor[(dx[V[n]]+Sum[V[s]V[n-s],{s,1,n-1}])/2],{n,1,M-1}];

Чтобы было удобнее смотреть, заменим u[n] на (Subscript[u, n]) и умножим на степень 2:

Do[Print[(Subscript[V, n]),"= ",(2^n)V[n]/.u[k_]:>(Subscript[u, k])],{n,1,M}]

Можно проверить, что построенный отрезок ряда действительно удовлетворяет уравнению (5), с точностью до членов (z^n) с n≥M:

v=-z+Sum[V[n]/(z^n),{n,1,M}];
Series[Expand[(v^2)+dx[v]],{z,Infinity,M-1}]

Прекрасно, и зачем это было нужно? А вот зачем: оказывается, что многочлены (Subscript[V, n]) являются плотностями законов сохранения для КдФ, то есть, существуют некоторые функции (Subscript[s, n]) от (Subscript[u, k]) такие, что

(7) 	(Subscript[D, t])((Subscript[V, n]))=(Subscript[D, x])((Subscript[s, n])).

Эти функции мы пока не знаем, но можем убедиться, что они есть, применив к этому равенству вариационную производную, тогда правая часть уйдёт в 0. Проверяем:

Table[
vard[star[V[n],u[3]+6u[0]u[1]]],
{n,1,M}]

Следует сделать одну оговорку. Если посмотреть внимательно на выписанную выше таблицу, можно заметить, что все (Subscript[V, n]) с четными номерами сами являются полными производными по x, ещё до дифференцирования по t:

Table[vard[V[n]],{n,2,M,2}]

Такие плотности не засчитываются, они отвечают тривиальным законам сохранения (Subscript[D, t])((Subscript[D, x])(a))=(Subscript[D, x])((Subscript[D, t])(a)), которые выполняются для любой функции a и любого эволюционного дифференцирования (Subscript[D, t]).  А с нечетными номерами все в порядке, это нетривиальные плотности, вариационная производная от них не есть тождественный 0:

Table[vard[V[n]],{n,1,M,2}]

Итак, можно сформулировать такие утверждения.

Теорема. 1) Все коэффициенты в разложении (6) являются плотностями законов сохранения для КдФ;
2) плотности с чётными номерами тривиальны; 3) плотности с нечётными номерами нетривиальны.

Доказательство 1) очень простое. Очевидно, что v это плотность закона сохранения для мКдФ. Так как преобразование Миуры выполняется тождественно по параметру λ=-(z^2), то все коэффициенты (Subscript[V, n]) также являются плотностями, но так как они выражены через (Subscript[u, n]), то это плотности уже для КдФ. Из одного-единственного, самого простого закона сохранения для одного уравнения мы вытащили целую бесконечную последовательность для другого. Кстати, для мКдФ мы тоже можем теперь получить бесконечную последовательность, делая прямую подстановку (4) в (Subscript[V, n]), но мы это пропустим.

Доказательство 2) и 3) оставляется в качестве упражнения (в файле с задачами есть подсказки).

Алгоритм интегрирования по частям

Алгоритм интегрирования по частям позволяет по заданной функции f от динамических переменных найти  функцию g такую, что f=(Subscript[D, x])(g) или убедиться, что такой функции не существует. Сейчас мы его опишем, запрограммируем и применим для восстановления функций (Subscript[s, n]) из (7).

Замечание. Слабым местом этого алгоритма, как и формулы E(f)=0, и вообще очень многих символьных вычислений является, казалось бы, безобидная операция сравнения с 0. Всегда ли легко проверить, что какое-то выражение равно или не равно 0? Для этого мы должны иметь эффективный алгоритм приведения этого выражения к некоторой однозначной канонической форме. Для полиномов это просто – нужно раскрыть скобки и собрать коэффициенты. Для рациональных функций уже труднее – нужно приводить к общему знаменателю. А если в выражение входят какие-то алгебраические функции, или функции, удовлетворяющие какому-то набору тождеств – задача может быть очень нетривиальной. Тем более, если в выражении десятки тысяч членов и общий объем исчисляется мегабайтами, что совсем не редкость. Это нужно иметь в виду, но сейчас у нас такие сложности не возникнут, так как мы работаем только с полиномами.

Цель алгоритма в том, чтобы найти представление заданной функции f в виде

f=(Subscript[D, x])(a) + b,

где либо b=0,  либо b=b(x,(Subscript[u, 0])) c (∂_0)(b)≠0,  либо b=b(x,(Subscript[u, 0]),...,(Subscript[u, k]))  c  (∂_k^2)(b)≠0.  Функция int[f] возвращает ответ в виде пары {a,b}. В ней мы уменьшаем порядок b, пока это возможно: если на каком-то шаге имеем b=p(x,(Subscript[u, 0]),...,(Subscript[u, k-1]))(Subscript[u, k])+q(x,(Subscript[u, 0]),...,(Subscript[u, k-1])) , то k можно уменьшить, занося под (Subscript[D, x]) интеграл от p по (Subscript[u, k-1]).

(* порядок по производным *)
ord[f_]:=Module[{v=vars[f]},If[v=={},-1,v[[-1]][[1]]]]

int[f_]:=Module[{a=0,b=Expand[f],c,k}, 
While[k=ord[b];     (* вычисляем k - порядок b *)
And[k>0,                   (* выполняем цикл, пока выполняется условие, что k>0 и *)
c=D[b,u[k]];     (* b линейна по старшей переменной *)
Expand[D[c,u[k]]]==0
],
c=Integrate[c,u[k-1]];   (* интегрируем по переменной, следующей по старшинству *)
a=Expand[a+c];                        (* пересчет a и b *)
b=Expand[b-dx[c]]
];
If[ord[b]==-1,    (* если в результате остаток не зависит от u, то это функция от х *)
{Expand[a+Integrate[b,x]],0},
{a,b}
]
]

Возьмем какое-нибудь выражение и посчитаем от него вариационную производную.

f=u[0]u[5]u[3]
vard[f]

Получилось не 0, значит f нельзя представить в виде полной производной. Применим интегрирование по частям. Получится пара, в которой второй элемент это то, что не лезет под знак (Subscript[D, x]), из-за нелинейности по старшей букве.

ab=int[f]
Expand[f-dx[ab[[1]]]-ab[[2]]]

А в этих примерах проинтегрировать получается.

f=g[u[0]]u[1]
vard[f]
int[f]

f=dx[(x^3)u[0]u[1]+u[2]u[7]u[4]]
vard[f]
int[f]

Проверяем (Subscript[V, n]) с чётными номерами. Да, действительно, в ответах вторая компонента равна 0, то есть всё занеслось под (Subscript[D, x]):

Table[int[(2^n)V[n]]/.u[k_]:>(Subscript[u, k]),{n,2,M,2}]

Проверяем (Subscript[V, n]) с нечётными номерами. Дифференцируем по t, в силу КдФ, и заносим под (Subscript[D, x]), опять вторая компонента равна 0, а первая это и есть искомые функции (Subscript[s, n]):

Table[int[star[(2^n)V[n],u[3]+6u[0]u[1]]]/.u[k_]:>(Subscript[u, k]),{n,2,M,2}]